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Abstract 

We consider surface-tension driven convection in a rotating fluid layer. For nearly 
insulating boundary conditions we derive a long-wave equation for the convection 
planform. Using a Galerkin method and direct numerical simulations we study the 
stability of the steady hexagonal patterns with respect to general side-band insta- 
bilities. In the presence of rotation steady and oscillatory instabilities are identified. 
One of them leads to stable, homogeneously oscillating hexagons. For sufficiently 
large rotation rates the stability balloon closes, rendering all steady hexagons un- 
stable and leading to spatio-temporal chaos. 
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In the study of the formation of patterns the origin of spatial and temporal dis- 
order is an important question. In recent years quite some progress has been made in 
the understanding of various states of spatio-temporal chaos that arise from stripe- 
like patterns, as they are found in the convection of fluids heated from below, or 
from spatially homogeneous oscillations. Much less is known about spatio-temporally 
chaotic states that are based on hexagonal patterns. We therefore investigate the 
classic surface-tension-driven Marangoni convection in a modification in which the 
system is rotated around a vertical axis. In buoyancy- driven roll convection such a 
rotation is known to induce chaotic dynamics through the Kiippers-Lortz instability. 
In order to obtain a systematic reduced description of the fluid system we consider 
the case of long-wave convection and derive the appropriate order-parameter equa- 
tion. Its stability analysis reveals regimes in which hexagons of all wavenumbers are 
unstable leading to a spatio-temporally chaotic state. In addition, homogeneously 
oscillating hexagon states are found. 



1 Introduction 

In the classic, buoyancy- driven Rayleigh-Benard convection rotation of the 
system has been found to lead to very interesting new dynamic patterns. The origin 
of the dynamics is the Kiippers-Lortz instability [1], which renders all straight roll 
patterns unstable with respect to a set of rolls rotated with respect to the original 
ones. It occurs above a critical rotation rate. In simple low-dimensional models the 
instability is found to a heteroclinic cycle connecting roll patterns rotated with re- 
spect to each other by 120° [2,3]. In experiments (e.g. [4,5]) as well as in numerical 
computations of the Navier-Stokes equations (e.g. [5]) and of model equations (e.g. 
[6]) more complex behavior is found. Patches of convection rolls of various orienta- 
tions arise and the roll orientation changes through the propagation of the fronts 
separating them or by the nucleation of rolls with a new orientation. 

The other classic convection system is surface-tension-driven Marangoni con- 
vection of a liquid layer in contact with another liquid or a gas layer (for a recent 
review see [7]). Two different instability mechanisms have been identified that drive 
this type of convection: a short-wave instability, which occurs independent of any 
deflections of the interface [8] and a long- wave instabihty, for which deformations of 
the interface are essential [9,10]. Not too far from threshold, the short-wave insta- 
bility typically leads to convection patterns with hexagonal planform (e.g. [11,12]), 
while the long- wave instability induces isolated dry spots or high spots [13,14]. The 
stability of the hexagonal short-wave pattern with respect to sideband instabilities 
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(the 'stability balloon') has been determined in weakly nonlinear [15] and in fully 
nonlinear numerical approaches [16]. Further above onset the short-wave hexagon 
pattern is often observed to undergo a transition to squares [17-19,11,20,12]. If the 
mean surface tension of the fluid is not too large the short-wave convection mode 
has to be coupled to a weakly damped long-wave mode that corresponds to the local 
layer height of the fluid [21,22]. A different type of long- wave instability has been 
identifled theoretically for the case when the fluid layer is bounded above and below 
by poor thermal conductors [8]. The mechanism for this instability is the same as 
for the short-wave instability. Due to the poor conductivity the critical wavenumber 
is, however, shifted to small values. Again one flnds that near onset the planform is 
usually hexagonal. It turns out that experimentally this limit of long wavelengths is 
difficult to achieve. Theoretically, however, it is a very useful hmit since it allows to 
derive simplified evolution equations systematically from the basic fiuid equations 
without the restriction to small amplitudes, although restrictions on the control 
parameter remain [23-25]. Moreover, compared to the weakly nonlinear Ginzburg- 
Landau equations the long-wave equations have the additional great advantage of 
preserving the isotropy of the system. This will play an important role in the present 
paper. 

With the interesting behavior ensuing from the Kiippers-Lortz instability, 
it is natural to ask whether rotation has also an interesting impact on hexagonal 
convection patterns. Within the general context of small- amplitude expansions and 
Ginzburg-Landau equations this question has been addressed previously [26-33]. 
These analyses were aimed at a general investigation of hexagonal patterns in rotat- 
ing systems and apply also to weakly non-Boussinesq Rayleigh-Benard convection. 
It was found that in a rotating system a new state of oscillating hexagons can ap- 
pear. This transition replaces then the transition from hexagons to rolls that is 
obtained in non-Boussinesq Rayleigh-Benard convection. In Marangoni convection 
this transition does not occur and no oscillating hexagons of this type are expected. 
Instead, in the absence of rotation a transition to squares is found [17-19]. The 
effect of rotation on this transition has not been studied yet. For the usual, steady 
hexagons the weakly nonlinear analysis indicates the possibility of persistent chaotic 
dynamics coexisting with the steady hexagons for the same parameter values [31]. 
However, these dynamics cannot be captured by the Ginzburg-Landau equations 
employed. To make reliable statements about this chaotic state, which is charac- 
terized by a spatially isotropic Fourier spectrum, equations have to be studied that 
preserve the isotropy of the system. A first step in this direction was the study of a 
Swift-Hohcnberg model, which was generalized to capture the breaking of the chiral 
symmetry that characterizes the rotation. A detailed analysis of the resulting equa- 
tion, which was not derived from any physical system, revealed that steady hexagons 
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of all wavenumbers can become unstable already for amplitudes below the transition 
to the oscillating hexagons, leaving a range in the control parameter over which no 
simple periodic solution is stable [34] and in which chaotic dynamics are found in 
numerical simulations. 

Motivated by the previous results of weakly nonlinear approaches to hexagons 
in rotating systems, our goal for the present paper is to obtain a systematic de- 
scription of Marangoni convection in a rotating system that is not limited to small 
amplitudes and that preserves the isotropy of the system. We therefore turn to the 
long-wave limit that is obtained for poorly conducting top and bottom boundaries 
and extend the derivation of [23,24] to the case of a rotating system. The poor 
heat conductivity is expressed in terms of a small Biot number Bi. While the true 
limit of long wavelengths {Bi <S 1) may not be attainable in experiments [14] we 
expect to gain valuable insight into the system for moderate wavelengths by ex- 
panding in Bi. The form of the resulting long-wave equation is the same as that 
found for buoyancy-driven convection with poorly conducting boundaries [35] and 
for buoyancy-driven convection in binary mixtures with positive separation ratio 
[36]. In the latter case the critical wavenumber goes to at a critical value of the 
separation ratio. In both cases quadratic terms arise only if the vertical symmetry is 
broken by imposing different boundary conditions at the top and the bottom plate 
(free-slip vs. no-slip or fixed temperature vs. fixed heat fiux). The focus of [35] is 
the effect of rotation on the system. In contrast to our work, however, there the 
competition between strictly periodic roll and square planforms is considered. The 
work is performed analytically in the weakly nonlinear regime, supplemented by a 
few simulations of a small system. 

In sec. 2 we introduce the basic fiuid equations and associated boundary con- 
ditions. The general linear stability analysis that identifies the long-wave regime is 
performed in sec. 3. The long-wave equation is derived in sec. 4. The stability and 
dynamics of the hexagon patterns is then studied numerically in sec. 5. The main 
results of the simulations are periodic and chaotic dynamics arising from secondary 
instabilities of the hexagons as well as from a resonant mode interaction. 



2 Formulation of the problem 
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2.1 Basic equations 



We consider a thin viscous and incompressible fluid layer heated from below 
with a free surface and rotating around a vertical axis. The mean surface tension is 

assTimcd to be large enough to suppress surface deformations even in the presence 
of rotation. The Navier-Stokes and the heat equation then read 



p{dt + Y ■ V) v = -Vp+ W|Qe^ X + 2pQe2 x v + rjAv - pge^, (1) 

V-v = 0, (2) 
(9t + v- V)r = KV2T. (3) 

In this paper we focus on fluid layers that are sufficiently thin that buoyancy can 
be ignored. The density is then taken to be constant, p = Pq. The dynamic viscosity 
f] is related to the kinematic viscosity i> through via i' = f]/po and k is the thermal 
diffusivity. The equations in dimensionless form can be obtained if the space, time, 
angular velocity, temperature, pressure and velocity fields are respectively divided 
by the constants d, (P/k,, k,/(P, AT, poKu/cP and n/d. Here d is the thickness of 
the fluid layer and AT is the temperature difference across the liquid layer in the 
conductive state. They are written as follows. 



1 V 

- (^t + u • V) u = -Vp + —V\n X r|2 + 2fi X u + V^u - Se^, (4) 

V-u = 0, (5) 
{dt + u-V)9 = V^9, (6) 

where V = f/nis the Prandtl number and Q = d^pg/{'qv) is the Galileo number. In 
this paper we consider the limit of infinite Prandtl number. 

We ehminate the pressure from the Navier-Stokes equations (4) by using the 
incompressibility condition (2). Taking the curl and the double curl of (1) we obtain 
[37] 



V^'y;-2Q9^C = 0, (7) 
V^C + 2J^9z«^ = 0, (8) 

where w is the ^-component of the velocity u = (-u, f , w) and is the 2;-component 
of the vorticity V x u. We also introduce the horizontal Laplacian Ah = {^xx'^^yy)- 
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Note that the rotation couples the vertical velocity and the vertical vorticity. Thus, 
while in the absence of rotation the vertical vorticity vanishes in many solutions of 
interest, with rotation the vertical vorticity is in general non-zero. 

The equations for the x and y-components of the velocity field are obtained 
considering that the velocity field u = (m, w) is expressed in terms of the poloidal 
and the toroidal stream functions (j) and as follows: u = {8x0^(1) + dy$,, dyd^cf) — 
dxC, —^h4>)- Then the equations for u and v are: 



- Ahu = d^d^w + dyC, (9) 
-Anv^dyd^w-dxC. (10) 

These equations together with (7)-(8) are equivalent to (4)-(5) [24]. 



2.2 Boundary conditions 



In order to obtain the velocity boundary condition at the top free surface 
we take into account thermocapillary effects and assume a linear dependence of the 
surface tension on the temperature, a — (Tq — ^{T — Tq). Furthermore, for sufficiently 
large mean surface tension and not too high rotation rates we can take the free 
surface to be flat. The boundary conditions at the upper free surface are then given 

by 



(dw du\ (dT\ 
(dw dv\ (dT\ 

After taking derivatives of (11) and (12) with respect to x and y, respectively, and 
addition of both expressions, one is led to 



(S^) = -7AhT. (13) 

Introducing the Marangoni number M = '^dAT / [pvK) the boundary condition can 
be written as 
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^ = MAh9 atz = l. (14) 

From the Marangoni equations (11-12) and from the underformabihty of the 
surface the other boundary conditions at the top are 



w = = 0, at^ = 1. (15) 
At the bottom we take no-shp boundary conditions for the velocity [38-40,35], 

w — ( — dzW = 0, at 2; = 0. (16) 

The critical wavenumber at which convection sets in depends sensitively on 
the temperature boundary conditions. Long wavelengths are obtained if the bound- 
aries are poor conductors. The Newton law at the bottom plate can be expressed 
as 

ka.T/b = hi{Tfb - Tb) at ^ = 0, (17) 

where T/b is the temperature of the fluid at the bottom plate and Th is the temper- 
ature of the plate at some distance below this interface. The thermal conductivity 
of the liquid is given by k and hi is the heat exchange coefficient at the bottom. 
We non-dimcnsionalize the temperature field by setting 9 = (T — Tfi)c)/AT, where 
Tfhc is the temperature of the fluid in the conductive state at the bottom plate and 
AT is again the temperature difference across the layer in the conductive state, i.e. 
AT = Tfbc — Tftc with Tftc the upper temperature of the fluid in the conductive 
state. In terms of 9 the boundary condition (17) is then given by 

dz9^Bii9-l atz = (18) 

with the Biot number Bii = hid/k. Here we have made use of the conservation of 
heat flux at the interface, hi{Tfi,c — T^) = — k/(i(T/bc — Tftc). In a similar way the 
temperature boundary condition at the top is found to be 

dz9 = -Bi29 -1-Bi2 at z = 1. (19) 

These expressions agree with those used in [41] when the limit of a single layer is 
taken (in their notation this limit corresponds to 0i = ©2, = 1 and a = 0). 



7 



3 Linear stability analysis 



In order to motivate the long-wavelength analysis we first discuss the linear 
stabihty of the conductive solution. The normal modes of the linearized system (6-8) 
together the boundary conditions (14-16), (18-19), take the form 

(u,^) = (u(^),^>))e''i-+^*, 

where A is the growth rate, q is the wavevector and r is the position vector (x, y). 
Since we are interested only in thermocapillary effects the Marangoni number is the 
sole control parameter. At the critical value Mc of the Marangoni number, deter- 
mined by a vanishing real part of A, the conductive state first becomes unstable. If 
the growth rate is A = at Mc the bifurcation is steady, while ii A — iuj for some 
nonzero real frequency cu, the bifurcation is oscillatory. In this work we examine only 
the steady bifurcation. It is known that for sufficiently large rotation rates an oscil- 
latory bifurcation may occur [42]. In the absence of rotation, the linearized system 
(6-8) yields the critical Marangoni number as a function of the wave number |q| 
and the Biot numbers Bii 2 



M — 8q ((— coshg -|- q coshg^ sinhgj (Bii + Bi2) 
-\- {^—qs\nh.q -\- coshg^ — coshgr^ BiiBi2 

-\-q^ coshg'^ — coshg — q^ sinhg^ /A^ (20) 
where, 

A = {coshq^ sinhg — sinhg — q^ coshg^ Bii -\- q^ coshg (21) 
— 2q^ sinhg + q coshg'^ — q coslig — g^ sinhg) . (22) 



As shown in Fig. la the critical wavenumber q^ vanishes for Bii = Bi2 = 0, 
while for Bii > and Bi2 > one has gc > 0. A series expansion of (20) around 
Bii = -6^2 = shows that gc = ((Sii -FSi2)405/27)^/^. Thus, for poorly conducting 
boundary conditions the wavenumber of the growing pattern is still small enough to 
allow the description by a long- wave model. In Fig. lb the marginal curve is shown 
for Bii = Bi2 = and two values of the rotation rate, Q = 1 and 1] = 5. We see that 
small values of the rotation rate do not modify the critical wavenumber. However for 
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M 




M 




Fig. 1. (a) Linear stability curve obtained from the full governing equations without rota- 
tion for Bii = Bi2 = (dashed line) and for Bii > and Bi2 > (solid line); (b) Linear 
stability analysis of the full governing equations for Bii = -6^2 = at = 1 (dashed line) 
and = 5 (solid line). 

> J7c = 3.7513 a minimum develops in the neutral curve even for Bii — Bi2 — 
and the critical wavenumber changes sharply from zero to a finite value. Thus, our 
study will be valid only for VL < Vt^. \n this regime the critical Marangoni number 
at g = as a function of VL is given by 



M,{n)\ 



q=0 



4 (2 cos{VQ) sin(v^) - sinh(2 V^)) 
4 sinh(v^) sin(v^) + cos(2 VU) - cosh(2 VU) ' 



(23) 



4 Derivation of the Long- Wave Equation 



Here we derive the long-wavelength model for rotating Marangoni convection 
and consider the system (3)-(10) together with the boundary conditions (14)-(15). 
Expanding 6 and u around the conductive solution, we rewrite the equations in 
terms of the deviations from the basic solution. To avoid excessive notation we keep 
the same notation 9 and u for deviations of the temperature and velocity fields and 
obtain 
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V^w -2nd,C^0, (24) 

V^C + 2nd,w = 0, (25) 
AHU + d^d,w + dyC = Q, (26) 

A„v + dyd,w~d,C = 0, (27) 
{dt + u-V)e = w + Ae. (28) 

The boundary conditions are: 

9;, - Bii9 = w = C = -u;^ = 0, at z = 0, (29) 
w^^ - MAhO = + 3120 = = = 0, at ^ = 1. (30) 

We define a small parameter e which measures distance from the critical Marangoni 
number by setting: 

M-Mc(l + eV). (31) 

Here is a parameter which allows to introduce changes in the Marangoni number 
independently of e. We assume that the heat loss is small, 

Bi = e%. (32) 

As discussed above, this implies that the dimensionless wavelength of the initial 
instability is long, of order e~^. This allow us to introduce the slow scales, 

X ^ex, y = ey. (33) 
The vertical coordinate z is not rescaled. The evolution occurs on the slow time T, 



t = e^T, (34) 

as suggested by the linear stability analysis [23]. These scalings of space and time 
suggest the scaling of the fields involved. Specifically, the Marangoni condition (30) 
suggests that the vertical velocity be taken to be small while the deviation of the 
temperature from the linear profile can be large. 



{w,e) = {e'^W,e). (35) 

Equations (26) and (27) point out the scaling for the x and y component of the 
velocity. 
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{u,v) = {eU,eV). 



(36) 



From the definition of tlie ;2-component of tlie vorticity and considering tlie previous 
scaling we have, 



C = e'Z. (37) 
With these definitions the equations (24)- (28) become: 

D^W + 2e^D^AHW + e^V% - 2fiDZ = 0, (38) 

D^Z + 2nDW + e^AHZ^O, (39) 

AhU + Wx, + Zy^O, (40) 

AhV + Wy.-Zx = Q, (41) 
e^Qr - D^e - A^^e - e^W = -e'^UQx - e^WQy - e^WDQ, (42) 

with D= d/dz. The boundary conditions are, 

W = DW = Z = D&-bie^e = 0, onz = 0, (43) 

W = B^W - Mc(l + eV) A^e = DZ = DG - 626^6 = 0, on ;2 = 1. (44) 



Proceeding as in Refs. [24,23], we expand the temperature as follows, 

e = eo(X, F, z, T) + ^Qi{X, F, ^, T) + 

and similarly for the other fields. Because in (38)- (44) only even powers of e appear, 
the expansion of the fields can also be restricted to even powers. If we insert these 
expansions in (38)-(44) we obtain at 0(e°) 



eo = F(X,F,T), (45) 

Wo = MAHFh{n,z), (46) 

Zo = MAHFf2in,z), (47) 

Uo = -M{FxBf2 + Fy^), (48) 

Vo = -MiFYBf2-Fxfi). (49) 



In contrast to the non-rotating case the vertical vorticity Z does not vanish and at 
each order we must solve two coupled inhomogeneous differential equations instead of 
only one as in [24]. In these expressions F{X,Y,T) is the planform function whose 
governing equation we are looking for. The functions /i and /2 are complicated 
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functions of Q and z. In the limit Q — > the function /2 vanishes and /i is the same 
as in [24] . 

At C(e^) we obtain, 

D^Oi = -AhF -Mfi AhF - MD/i \VhF\^. 

The critical Marangoni number at A; = is obtained from the solvability condition 
in agreement with the expression (23) obtained in the linear stability analysis. For 
the fields we obtain 

Bi = MQ, z)AhF + f,{Q, z)\VhF\^ + F^iX, Y, T), (50) 
W^ = h{n, z)AhF + fe{n, z)A%F + fr{n, z)Ah\VhF\^ + 

fs{n,z)AHF2, (51) 
Z, = fg{Q, z)AhF + z)A]jF + /n(l^, z)Ah\VhF\' + 

/i2(l^,^)A^F2, (52) 



Jl2[}L,Z)L\H^2, 

Ui — —fg Fy — fio AffFy — fn (| V_f/Fp)y — /12 F2^y — f^^z Fx — 
h,z AhFx — f7,z{\^ hF\^)x — f8,zF2,X, 

Vi = /g Fx + /lo AhFx + fu (I Vj^FHx + /12 F2,x - h, " 

h,z AnFy — /7,z(|VifF|^)Y — f8,zF2,Y- 



(53) 
(54) 



Here F2 is the 0(e^)-correction to the planform function F. The functions f^, f^, 
fbi /g) /t) /s) /q) /lO) /ii si'iid /12 are again complicated functions of z and Vt. The 
functions /a and are related to integrals of /i with respect to z. It turns out 
that after the integration these functions contain terms that are singular at Q = 
and do not depend on z in this limit. In order to cancel these singularities and to 
recover the expressions of [24] in the limit — > 0, the planform function F2, which 
is undetermined at this order, must contain the same singular constants with the 
opposite sign. Since the integral expressions are too complicate to present here, this 
is best illustrated by the simple example 



j iim.{Q,z)dz = cos{9,z) + C. (55) 

This integral is not finite in the limit ^ 0. If a finite value is required for it, the 
arbitrary constant C must be chosen as l/Vt. Analogously, in and constants 
appear that need to be chosen appropriately. Since they do not depend on z they 
can be absorbed in F2 and therefore do not affect the long-wave equation for F. 
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Fig. 2. The constants A, ci, C2, C3, C4 and C5 plotted as functions of the rotation rate Q,. 

Finally, at 0{e'^) the solvability condition for 0^62 yields the required equa- 
tion for the planform function F{X, Y, T), 



Ft = -(61 + b2)F - i^AhF - XAj^F - ciV h ■ {VhF^VhF - 

C2VH ■ AhFVhF - czAh\VhF\^ - c^e, ■ VhF x Vh{AF} 
+C5e,-VjjFx Vj^{(VifF)2}. 



(56) 



Here fi is the control parameter (cf. (31)) The constants A, Ci, C2, C3, C4 and C5 
" — complicated functions of the rotation rate ^2. We plot them in Fig. 2. A series 



are 
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expansion of these functions around = 0, shows that in the hmit Q — > 0, the same 
equation as in [24] is obtained: 



5 Numerical results 

We now turn to the stabihty and dynamics of hexagon patterns described by 
(56). One of our goals is to identify parameter regimes in which regular hexagon pat- 
terns become unstable at all wavenumbers. In the weakly nonlinear case hexagons 
can be described by three coupled Ginzburg-Landau equations. Within this frame- 
work their stability and dynamics has been studied in detail without and with ro- 
tation [43-46,31,32,47]. It has been found that for small amplitudes there is always 
a range of wavenumbers over which the hexagons are stable. In a simple Swift- 
Hohcnberg model, situations were identified in which somewhat further above onset 
the stability balloon closes, leading to persistent chaotic dynamics [34]. In this pa- 
per we therefore do not reduce (56) to small-amplitude Ginzburg-Landau equations. 
Instead determine the stability of the hexagon solutions of (56) and the resulting 
dynamics numerically without any restriction to small amplitudes. 

The stability of the hexagon patterns is determined with a Fourier-Galcrkin 
code on a hexagonal lattice. We then determine their stability with respect to 
side-band perturbations using a Floquet analysis with a Floquct parameter p = 
p(cos^, sin^). Due to the isotropy of the system and the discrete rotational sym- 
metry of the hexagons we can restrict our analysis to perturbation wavevectors 
subtending an angle 9 of up to ±30° relative to the hexagon modes. It turns out 
that accurate solutions of (56) require a large number of Fourier modes. One measure 
of the accuracy of a solution is obtained with the stability analysis by determining 
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n 


q 


n = 7 


n= 11 


n = 15 


0.07 


3.0 


0.98 


-0.16 • 10-^ 


-0.28 • 10-5 


-0.56 • 10-'^ 


0.17 


2.5 


0.95 


-0.25 • 10"^ 


-0.53 • 10-2 


-0.60 • 10-3 



Table 1 

Translation eigenvalue cr(p = 0) as a function of the truncation parameter n of the Fourier 
expansion for 6i + 62 = 0.025 illustrating the poor convergence for larger amplitudes. 




I 0g I I I I I I I I I I I I I I I'T' I I I I I I I I I I I I I I I I 

* 0.65 0.70 0.75 0.80 0.85 0.90 0.95 



Wavenumber q 



Fig. 3. Stability limits of the steady hexagons for = and 61 + 62 = 0.025. The dashed 
line gives the neutral curve. The stability limit is given by a steady long-wave instability 
(solid line). The plusses and crosses indicate stable and unstable solutions according to 
the numerical simulation. 



the magnitude of the eigenvalue corresponding to the translation mode given by the 
Floquet parameter p = 0. For the exact solution one has o'(p) = 0, while numerically 
truncated solutions will have o'(p) 7^ 0. The dependence of o"(p) as a function of the 
number n of harmonics retained is illustrated in table 1, where n = 1 corresponds 
to the basic hexagon solution. 

Fig. 3 shows the resulting stabihty hmits in the absence of rotation (Q = 0) 
for 61 + 62 = 0.025. They are given by steady long-wave instabilities with p ^ 
1 (cf. [43,44,47]). Numerical simulations for selected parameter values (marked by 
plusses and crosses in fig.3) confirm that the instabilities lead to the formation of 
defects which then change the wavevector of the pattern. 
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Fig. 4. Stability limits of the steady hexagons for = 3 and 6i + 62 = 0.025. Hexagons 
are only stable inside the closed loop. The neutral and saddle-node lines are given by thin 
dashed and solid lines. The symbols identify the character of the relevant instabilities: 
long- wave steady (open circles), short-wave steady (open left-triangles), long- wave oscil- 
latory (solid circles), homogeneously oscillatory (solid squares). The plusses and crosses 
denote simulations in which the steady hexagons were stable/unstable (see text). 

Rotation strongly modifies the stability limits as shown in fig.4 for Q = 3 
and bi + b2 — 0.025. The most striking aspect is that the stabihty balloon closes 
for II > 0.07. This closing is due to a steady short-wave instability (triangles) which 
supersedes the usual steady long-wave instability (open circles) for larger ^. In the 
short-wave regime the largest growth rates are obtained for p 0.6q and 9 ~ 
25° — 30°. The appearance of a short-wave instability is somewhat similar to the 
result found in the weakly nonlinear analysis. There it has been found that in the 
presence of rotation a steady or oscillatory instability arises at finite modulation 
wavenumbers and preempts the long-wave oscillatory instability that emerges from 
the interaction of the two long- wave phase modes [31]. 

A further interesting feature of the stability limits for Q = 3 is the oscillatory 
instability that hmits the stability balloon on the low-q side (solid symbols). For low 
values of /i it is a long-wave side-band instability, i.e. its growth rate and frequency 
vanish for p = (solid circles). For larger /i, however, the largest growth rate occurs 
at p = and the hexagon pattern is unstable to a homogeneous oscillatory mode 
(solid squares). In the weakly nonlinear analysis such a homogeneous oscillatory 
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Fig. 5. Maximal value of the control parameter for which there are still steady hexagons 
(solid line) and onset of hexagons (dashed line) (61 + 62 = 0.025). 

instability has only been identified in the context of the transition from hexagons 
to rolls [26,27,31]. No such connection is apparent in the present system. 

Anticipating interesting dynamics for values of /i above the closing of the 
stability balloon, we have studied its dependence on the rotation rate. As a function 
of Q, fig. 5 shows the maximal control parameter fi for which there are still stable 
steady hexagons (solid squares). To indicate the /i-range of stable hexagons the 
critical control parameter is also given (dashed line). With decreasing Q that 
/i-range grows rapidly. To resolve the stability limits for yet smaller values of f2 a 
very large number of modes is needed; = 19 is not sufficient (cf. table 1). 

The nonlinear evolution ensuing from the linear instabilities is investigated by 
simulating (56) numerically with a pseudo-spectral Fourier code using an integrating- 
factor, fourth-order Runge-Kutta scheme. It turns out that for larger values of /i the 
nonlinear gradient terms in (56) can under certain conditions lead to an apparently 
spurious blow-up of the numerical solution in which very high wavenumbers are 
excited [48][^ However, if a sufficiently smooth initial condition was used no such 

^ Even though the numerical time-stepping routine may lead to a blow-up, no change 
in the behavior of the stationary states is found with the Galerkin code for those pa- 
rameters. Apparently, in these cases the Fourier representation of the derivatives in the 
time-stepping code does not provide sufficient numerical dissipation. By introducing an 
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Fig. 6. Snapshots of oscillating hexagons for 17 = 3 and 6i + 62 = 0.025, n = 0.060, and 
q = 0.93 at times t = 0, t = ^T, t = OAST, and t = fT, where the period T = 140 (time 
progressing left to right, top to bottom; red denotes positive, blue negative values). 

numerical divergencies occurred. As already indicated earlier (cf. table 1) the proper 
resolution of the patterns requires a large number of Fourier modes. In all simulations 
we have used a rectangular computational domain of size Lx = 2L and Ly = 2L/V3 
to accommodate undeformed hexagon patterns. Specifically, for L = N2TT/q this 
domain allows N periods of the pattern in the x- and y-direction. 

To test the stability results we have performed simulations near the stability 
limits depicted in figs. 3,4. The plusses denote parameters for which a small modula- 
tion of the pattern did not destabilize it, while at the crosses the hexagons became 
unstable. These simulations were performed in systems of size L = 4 ■ 27r/g. For this 
aspect ratio some of the instabilities are still slightly suppressed shifting the stabil- 
ity limits somewhat. The limits obtained in the simulations are, however, consistent 
with the stability results if the modulation wavenumbers p are restricted to those 
compatible with the computational domain used in the simulations. 

We first address the nonlinear evolution of the homogeneous oscillatory insta- 
bility (solid squares). Simulations for /i = 0.060 and q = 0.93 show that this stability 
limit is due to a supercritical Hopf bifurcation, i.e. the spatially homogeneous os- 
cillations saturate and lead to a state of oscillating hexagons. They are similar in 
appearance to the oscillating hexagons that arise in the Hopf bifurcation that re- 
places the instability of hexagons to the mixed mode, which is associated with the 
transition to rolls. As mentioned earlier, this transition to the mixed mode and to 
rolls has not been seen in the present system. A time sequence of snapshots of the 
oscillating hexagons is shown in fig.fQ. 

Now we turn to the nonlinear evolution arising from the steady short-wave 
instability (open left-triangles in fig. 4). In fig. 7 and fig. 8 a sequence of the pattern 
and the temporal evolution of two modes is shown, respectively in a system of size 

artificial dissipation term these divergencies can be avoided. To obtain sufficiently smooth 
starting solutions we have at times used a filter that affects only wavenumbers above 
a certain high wavenumber and corresponds then to a linear 8th-order derivative term. 
However, all the simulations presented in this paper have been performed with that filter 
turned off. 

^ A movie of the temporal evolution of the oscillating hexagons can be found at 
http : //www. esam. northwestern. edu/riecke/research/Marangoni/marangoni .html . 
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see fig7artc.jpg 

Fig. 7. Snapshots of typical patterns obtained for $7 = 3 and 61 + 62 = 0.025 at 
Lx = 20.8,// = 0.07 at times 6299, 6349, 6449, 6499, 6549, 6599, 6649, 6749 (left to 
right and top to bottom) 
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Fig. 8. Time series obtained for = 3 and 61 + 62 = 0.025 at Lx = 20.8, /x = 0.07. It shows 
the energy S associated with the Fourier modes with wavevectors (27r/20.8, 2 • 27r\/3/20.8) 
with thin and thick solid lines, respectively. Note the long onset time for the instability. 

fig9artc.jpg 

Fig. 9. Typical patterns obtained for = 3 and 61 + 62 = 0.025 at Lx = 48, ^ = 0.08 at 
times 2699, 3749, 3899, 11549 (left to right and top to bottom). 

Lx = 20.8 for = 0.07 (using 128 x 128 modes). For this value of the control 
parameter no stable hexagon patterns exist in the infinite system. In this case the 
temporal evolution does not become periodic for the time period investigated. 

In fig. 9 typical snapshots for the evolution of the pattern in a system of 
size Lx = 48 are shown for = 0.08 (using 256 x 256 modes). The associated 
time series for two Fourier modes are shown in fig. 10. Again, for these parameter 
values no steady stable hexagons exist in large systems resulting in persistent chaotic 
dynamics.[3 While the pattern is spatially quite disordered a tendency to form locally 
square arrangements can be observed. This presumably reflects the fact that in the 
absence of the quadratic terms in (56) even in the rotating system square patterns 
are preferred over roll patterns [35]. 

Finally, in Figs. 11, 12 we show the temporal evolution in a quite small system 

^ A movie of the chaotic dynamics for fi = 0.075 and Lx = 48 can be found at 
http : //www. esam. northwestern. edu/riecke/research/Marangoni/marangoni .html. 
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see figl0.jpg 

Fig. 10. Time scries obtained for f7 = 3 and 61 + 62 = 0.025 at Lx = 48, = 0.08. 
The energy S associated with the amplitude of the Fourier modes with wavevector 
(2 • 27r/48,4 • 27rV3/48) and (0,4 • 27rV3/48) are shown with thin and thick solid lines, 
respectively. 

see figllartc.jpg 

Fig. 11. Patterns in a cycle obtained for = 3 and 5i + 62 = 0.025 at Lx = 17, /x = 0.058 
at times 4079.2, 4159.2, 4239.2, 4319.2, 4399.2, 4479.2, 4559.2, 4639.2 (left to right and 
top to bottom). 




2000 4000 6000 8000 10000 ^ 12000 



Fig. 12. Time series obtained for J] = 3 and 61 + 62 = 0.025 at Lx = 17, /x = 0.058. It 
shows the energy S associated with the Fourier mode with wavevector (27r/17, 2-27r\/3/17). 

of size Lx = 17 for fi = 0.058 (using 64 x 64 Fourier modes). For such a small system 
the stability limits as shown in figs. 4 do not really apply. The initial condition was 
a regular hexagon pattern with wavenumber q — 1.478. After a transient during 
which defects appeared and changed the number of cells the time dependence of the 
pattern becomes quite regular. The time scries (cf. fig. 12) is, however, not strictly 
periodic but shows small 1%) irregular variations of the oscillation amplitude. A 
closer inspection of the neutral curve, which in terms of the system size and of the 
number of wavelengths contained in the system is given by 




shows that the parameters for this simulation are close to the codimension-two 
point at which hexagons with q — 0.74 and hexagons with q — 1.48 first arise. This 
is illustrated in fig. 13 which gives the neutral curve (57) for various values of A^. In 
this small system we associate the origin of the dynamics with the low-dimensional 
interaction of the two modes associated with the codimension-two point. 
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Fig. 13. The neutral curve ^{Lx) plotted for different values of N, the number of wave- 
lengths contained in the system. The parameters are 61 + 62 = 0.025 and A = 0.0209 
(which corresponds to a rotation rate of = 3). Squares stand for numerical simulations 
performed in different boxes at different /x values. 



6 Conclusion 



In this paper we have considered Marangoni convection in the presence of 
rotation. In order to derive systematically a set of reduced equations that - unlike 
the Ginzburg-Landau equations - preserve the isotropy of the system we have con- 
sidered poorly conducting boundary conditions. Under these conditions convection 
sets in first at long wavelengths, which allows the derivation of a long-wave equa- 
tion. Extending previous analyses [24,49], we have derived the long-wave equation 
in the presence of rotation and obtained the coefficients as a function of the rotation 
rate Vl. Using a Galerkin approach we have then determined the linear stability of 
steady hexagonal convection with respect to general side-band perturbations without 
and with rotation. We have employed direct numerical simulations of the long-wave 
equation to investigate the fate of these instabilities. 

In the absence of rotation the band of stable steady hexagons in large systems 
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is determined solely by a steady long-wave instability for the range of the control 
parameter considered. With rotation, however, steady long- and short-wave instabil- 
ities as well as oscillatory instabilities arise. The latter can be a long-wave side-band 
instability or an instability that is spatially homogeneous at its onset. Interestingly, 
the latter is supercritical and leads to stable, spatially homogeneous oscillations 
of the hexagons. This state resembles the oscillating-hexagon branch that arises in 
non-Boussinesq buoyancy-driven convection in connection with the transition from 
hexagon to roll convection [26,27,31]. In the present system no transition to rolls 
has been identified and in experiments on Marangoni convection a transition from 
hexagons to squares is observed instead [17,19]. We find temporally almost periodic 
and spatially more complex states that appear to arise from a resonant q : 2q mode 
interaction. This and other resonances could be studied in future work. Most inter- 
esting is the fact that for sufficiently large rotation rates the stability limits close for 
larger values of the control parameter. In the numerical simulations this is shown to 
lead to spatially and temporally chaotic patterns. 

So far, there have been very few experimental studies of spatio-temporal chaos 
arising from a hexagonal planform [50]. Therefore, it would be particularly inter- 
esting to see whether the complex dynamics of the disordered hexagonal patterns 
predicted by the long-wave model (56) can be obtained in experiments. It is yet to 
be seen to what extent in the current system the deformation of the free surface by 
the centrifugal force modifies the results presented here. Similarly, the derivation of 
the long-wave equation depends on the poor thermal conductivity of the horizontal 
boundaries. For experimentally realizable boundary conditions the wavelength of the 
convecting state may not be as large [14]. As is often the case, however, the insta- 
bilities predicted by the long-wave equation (56) may persist qualitatively beyond 
the asymptotic long- wave limit. A study of this question will presumably require 
the derivation of the Ginzburg-Landau equations with nonlinear gradient terms or, 
more hkely due to the finite amphtudes at the closing of the stabihty balloon, a 
numerical analysis of the full fiuid equations (cf.[16]). 
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